Quantum Phase Transition, 0(3) Universality Class and Phase Diagram of Spin- 1/2 Heisenberg 
Antiferromagnet on Distorted Honeycomb Lattice: A Tensor Renormalization Group Study 
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The spin-1/2 Heisenberg antiferromagnet on the distorted honeycomb (DHC) lattice is studied by means 
of the tensor renormahzation group method. It is unveiled that the system has a quantum phase transition 
of second-order between the gapped quantum dimer phase and a collinear Neel phase at the critical point of 
coupling ratio Oc — 0.54, where the quantum critical exponents v ~ 0.69(2) and 7 ~ 1.363(8) are obtained. 
The quantum criticality is found to fall into the 0(3) universality class. A ground-state phase diagram in the 
field-coupling ratio plane is proposed, where the phases such as the dimer, semi-classical Neel, and polarized 
phases are identified. A link between the present spin system to the boson Hubbard model on the DHC lattice 
is also discussed. 

PACS numbers: 75.10.Jm, 75.40.Mg, 64.70.Tg, 05.30.Jp 



I. INTRODUCTION 

Since the discovery of high temperature superconductors, 
the two-dimensional (2D) Heisenberg models have received 
particular attention in the past decades. Several numerical 
works show that the spin-1/2 isotropic Heisenberg antifer- 
romagnet (HAF) on a square lattice exhibits an AF long- 
range order (LRO) in the ground state^'^, although a math- 
ematically rigorous proof still lacks now. Various methods 
(e.g. the spin wave analysis"^, different numerical techniques^, 
etc.) were also applied to investigate the properties of this 
model. When the bond anisotropy is introduced, the mag- 
netic order-disorder quantum phase transition (QPT) can be 
identified^—. Another intriguing 2D bipartite lattice — the 
honeycomb (HC) lattice has also been studied with vari- 
ous methods, such as quantum Monte Carlo (QMC)^, series 
expansion'', spin wave'^ and newly proposed tensor renor- 
malization group (TRG) method'^. These investigations show 
that owing to the lowest coordinates among 2D lattices, the 
system is more affected by quantum fluctuations, giving rise 
to the spontaneous magnetization per site of the spin-1/2 HAF 
model on this lattice smaller than that on a square lattice. 

Recently, people have obtained a number of magnetic ma- 
terials with distorted honeycomb (DHC) lattices, such as 
MnPSs and FePS3'^ CU2/3V1/3OA NaaCuaSbOe'^ and 
Mn[CioH6(OH)(COO)]2X 2H20ii, where the magnetic ions 
(e.g. Cu+^ and Mn+^) form an HC lattice with different 
nearest-neighbor (NN) bonds. The magnetic properties of 
these materials have been investigated experimentally. How- 
ever, the theoretical studies on the HAF model on the DHC 
lattice are still sparse. There is a recent work that explores the 
ground state properties of the Heisenberg model on a DHC 
lattice by QMC calculations"^, and the order-disorder QPT in 
zero magnetic field has been observed. Nevertheless, the mag- 
netic properties of the model on such a DHC lattice in nonzero 
external fields are not yet seen in literature. Therefore, in or- 
der to understand the experimental observations profoundly, 
it should pay more theoretical attention on the HAF model on 
the DHC lattice. 

In this paper, by means of the newly developed TRG 
method and a variational analysis, we shall study the ground 



FIG. 1: (Color online) The distorted honeycomb lattice can be 
viewed as a tensor network, where Ta (Tb) are the tensors located 
on A (B) sublattice (indicated by different symbols), and each tensor 
has three bond indices labeled by x, y, and z. 



state properties of the spin-1/2 HAF model on the DHC lat- 
tice with and without magnetic fields. On one hand, our pri- 
mary purpose is to understand the ground state properties of 
the model under consideration in the presence of a magnetic 
field, and on the other hand, we also want to test the accuracy 
of TRG methods in a more extensive range rather than a spa- 
tially bond isotropic case by comparing our calculated results 
with the previous studies on HC and DHC lattices using other 
methods. Owing to the limitation of the present TRG scheme, 
a finite-temperature calculation on the Heisenberg model is 
still not feasible now, which makes it impossible to compare 
directly the calculated results with the experiments. Our study 
shows that the spin-1/2 Heisenberg model on a DHC lattice 
has a second-order QPT with respect to the bond coupling ra- 
tio, that is determined to fall into the 0(3) universality class 
by identifying the two quantum critical exponents. A phase di- 
agram separating the dimer, polarized, canted Neel and Neel 
phases is also proposed. The magnetic properties in the pres- 
ence of a magnetic field are obtained, where some interesting 
behaviors are observed. In addition, the TRG method has been 
verified through this present model to give a good agreement 
with most of the previous studies, but it might overestimate 
slightly the spontaneous sublattice magnetization per site on 
the HC lattice, implying that the TRG method may still need 
more works for improvements. 

The other parts of this paper is organized as follows. In Sec. 



II, the model Hamiltonian and the TRG method will be intro- 
duced. In Sec. Ill, the magnetic properties of the model under 
interest in the absence of a magnetic field will be reexamined 
with the TRG algorithm. In Sec. IV, the ground state proper- 
ties in the presence of a magnetic field will be explored, and a 
phase diagram will be proposed. In Sec. V, by invoking a bo- 
son mapping, the present model is mapped onto the 2D boson 
Hubbard model, whose properties will be briefly discussed. 
Finally, a conclusion will be given. 



II. MODEL HAMILTONIAN AND NUMERICAL METHOD 
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The system under interest is schematically depicted in Fig. 
1 . The Hamiltonian reads 



H = J 
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(1) 

where S^ denotes the spin-1/2 operator at site i, {i,j)x labels 
the NN spins along the rungs {x bond), {i,j)yz means the NN 
spins along the zigzag directions (y and 2 bonds), J is the 
interaction on the x bond, J' is the coupling on the y or z 
bond, h and hs stand for the uniform and staggered magnetic 
fields, respectively, and e^ = +1 when i on A sublattice and 
— 1 on i? sublattice. We introduce for convenience a bond 
coupling ratio a = J' /J, and take J as an energy scale. 

To explore the ground state properties of the present system 
we shall utilize the TRG method. This numerical algorithm 
was first introduced to calculate the thermodynamic proper- 
ties of the 2D classical models'^, and then generalized to ob- 
tain the expectation values of observables in the quantum state 
with the tensor product wave functions (e.g. RefUSl) on bipar- 
tite lattices given by 
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where Ta (Tb) represents the tensor located on A (B) sub- 
lattice, over which the indices i and j run, and the summa- 
tion over all bond indices x, y, z is from 1 to the bond di- 
mension D. According to the TRG algorithm, the ground 
state wave function and energy can be directly obtained by 
using trial wave functions of tensor product form^". This 
variational scheme, however, is not so efficient that makes 
the achievable bond dimension D not larger than 3 in gen- 
eral due to huge variational parameter space. Recently, it was 
improved by combining the infinite time-evolving block dec- 
imation (iTEBD)^' and TRG method to determine the ground 
state and to get the expectation values of local observables'^. 
This alternative algorithm appears to be accurate and efficient, 
in which the available bond dimension D can reach as large as 
8, and the calculated results agree well with those obtained by 
other methods. It has been applied to study the spin flop tran- 
sition of spin-1/2 XXZ model on a square lattice^^. In the 
following, we shall adopt this novel scheme to calculate the 
physical quantities of the spin-1/2 HAF system on the DHC 
lattice. 



FIG. 2: (Color online) (a) The spontaneous sublattice magnetiza- 
tion per site \'m^jy(g-^ \ as a function of the bond ratio a. for different 
bond dimensions D = 3,4, 5, 6. The dashed dot line represents 
the mean-field result [Eq. ([5]l]. The inset gives the derivative of 
j"^z,A(s)! as a function of a, where a discontinuity at a^ ~ 0.54 is 
observed, (b) The ground state energy per site e and its derivatives 
(insets: first- and second-order) versus q, where the discontinuity in 
d^e/do? against a indicates a QPT of second-order 
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FIG. 3: (Color online) The magnetization per site as a function of 
magnetic field h. The insets show the susceptibility x (upper panel) 
and the transverse component of sublattice magnetization per site 
\''T>,x.A{B) I (lower panel) as functions of h. (a) a = 0.3; (b) a = 0.6, 
where D = 5 and /is = for both. The fitting curves in (a) shows a 
nearly linear behavior of m^, with /? ~ 0.98(1). The number in the 
parenthesis denotes numerical fitting error hereafter 



In our practical calculations, during the iterative projections 
by evolution operator ie^'^^) along the imaginary time r axis, 
we first start with a step 5t — 10^'^, and then diminish it 
gradually to 5t = 10^^. The total number of iterations is 
taken as about 10^ ^ 10®, where D ~ 5 or 6 is generally 
chosen. The convergence is always checked, as shown in Fig. 
12 a) for different bond dimensions D. 



III. MAGNETIC ORDER-DISORDER TRANSITION 

Let us first consider the case in absence of a magnetic 
field (/i = and hg — 0). When a = 0, the spins are 
coupled only by J along x bonds, and the ground state is 
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per bond, which is usually termed as a dimer state. This disor- 
dered ground state is protected by a finite spin singlet-triplet 
gap. When J' is set in, but a is still small, one may conceive 
that the system may retain in the dimer phase^. This is con- 
firmed in Fig. |2{a), where the spontaneous sublattice magneti- 
zation per site, m^^Ais) = jrj^ J2^eA{B)iS^) with Na(b) 
the total number of sublattice A(B) sites, as a function of a 
is presented for Z) = 3, 4, 5, 6. It is seen that there exists a 
critical ratio ac below which m^^AiB) vanishes, showing that 
the ground state for a < ac is disordered and dominated by 
quantum fluctuations. For a > ac, the ground state of the 
system has an AFLRO owing to a spontaneous SU{2) sym- 
metry breaking. To determine the value of ac with accuracy, 
we have calculated the derivative of Im^AiB) I with respect to 
a for D = 6, and found a discontinuity at a — ac, where ac 
can be readily determined as 0.54. 

To confirm if the transition occurring at ac is a QPT, we 
have also studied the ground state energy per site e as a func- 
tion of a. The results are presented in Fig. lUb). It may be 
observed that with increasing a, both e and its first derivative 
de/da versus a decrease continuously, but the second deriva- 
tive d^e/da^ against a shows a discontinuity at ac — 0.54, 
as seen from the inset of Fig. |2jb). This feature characterizes 
a typical QPT of second-order 

It should be noted that, as shown in Fig. |3a), we get 
1*^2. A(B) I — 0.32 (D — 6) for a — 1, which appears to over- 
estimate the spontaneous sublattice magnetization per site on 
the HC lattice in comparison to the recent stochastic series 
expansion QMC resuh 0.2681(8) (Ref. IH), flie "world fine" 
QMC result 0.22 (Ref. ^, series expansion 0.27 (Ref. [01), 
and the spin-wave result 0.24 (Ref. 12). Such a discrepancy 
on mz^A{B) has also been noted in Ref. |2J, where they re- 
ported the sublattice magnetization per site on the HC lattice 
to be 0.3098 for D = 16 by the TRG calculations. This slight 
discrepancy on the sublattice magnetization per site may come 
from the underestimation of quantum fluctuations in the ab- 
sence of a magnetic field in the assumption of the tensor prod- 
uct state employed in the TRG method, because the isotropic 
system at a = 1 is gapless and has long-range correlations, 
where the quantum fluctuations may be strong. However, the 
ground state energy per site we obtained (-0.5465) for a — 1 
(D — 6) is quite consistent with those obtained for the HC lat- 
tice by other methods, e.g. the QMC result -0.5450 (Ref. 0), 
series expansion -0.5443 (Ref. II), and the spin-wave result 
-0.5489 (Ref. [l2l) . Therefore, as the critical point determined 
by the spontaneous sublattice magnetization [Fig. |2](a)] co- 
incides with that obtained from the singularity of the ground 
state energy [Fig. |2](b)], it shows that the critical point is de- 
termined with rather assurance. In the presence of a magnetic 
field, since the quantum fluctuations are much suppressed, the 
results given by the TRG method should be reliable. This 
can also be validated in the following studies, including the 
attained linear behavior of the magnetization curves immedi- 
ately above the critical magnetic field he [Fig. [3] (a)] and the 
verified relation Hsat = 2S'(2a + 1) of the saturation line that 
separates the canted Neel phase and polarized phase in Fig. 
m which is also consistent with that derived from the classi- 
cal energy of spin wave analysis by a variational schemei^. 



Another fact is that when the critical exponents of i/ and 7 are 
determined (in Sec. IV), the critical point is approached from 
the dimer phase, and hence it is independent of the magnitude 
of spontaneous sublattice magnetization in the Neel phase. 
The obtained v and 7 agree well with previous calculations 
and theoretical predictions, showing again that TRG method 
is fairly feasible for the present case. 

In order to examine our numerical results, we perform a 
mean-field treatment in terms of a simple variational trial 
wave function 
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that was applied to describe both disordered and ordered 
phases on a dimerized square lattica^i and in a bilayer 
system^^, where the lattice site i and i + a; are NN sites along 
X bonds, and i is a variational parameter and interpolates be- 
tween a singlet collection (t ~ 1) and a classical Neel state 
{t = 0). Substituting Eq. ([I} (with h, hs = 0) and Eq. OJ into 
Eyar = {'^var\H\'i'yar) and minimizing it with respect to t, 
we obtain an upper bound for the ground state energy as 



-3/8, 



for a < 0.5 



' ^ l/a + 4a + 2, for a > 0.5 (4) 

16 

where N — Na + Nb is the total number of lattice sites. Ob- 
viously, Eyar is singular at ac,var — 0.5, showing ac,var, 
that is close to ac = 0.54, may be a transition point. 
The variational sublattice magnetization per site niyar = 
{1/Na) J2ieA{'^var\Sf\^yar) Can be obtained by 

f 0, for a < 0.5 

"^™''~ j -^/-l/a^^+4:, for a > 0.5 (5) 

which indicates that the derivative of niyar is discontinuous 
at a = 0.5, suggesting a discorder-order phase transition at 
ctc.var- A comparison of rriyar to the TRG result is given in 
Fig. |2ja), where rriyar shows a behavior similar to the TRG 
results. 



IV. GROUND-STATE PHASE DIAGRAM AND CRITICAL 
EXPONENTS IN PRESENCE OF A MAGNETIC FIELD 

Now we turn the uniform magnetic field h on. The magne- 
tization curves for different a are given in Fig. [3] It is clear 
that for a < ac, there exists a magnetization plateau with 
TTiz = 0, where m^ = (l/N) J^ii^i) i^ the magnetization 
per site, implying the existence of a finite spin gap, as shown 
in Fig. Oa) for a = 0.3. At a = 0, such a gap is nothing but 
the spin singlet-triplet gap which equals J. For < a < ac, 
the spin gap will decrease and vanish eventually at a = ac- 
For a given a < ac, the spin gap closes at a critical field he- 
For h > he, we find that the magnetization depends almost 
linearly on the magnetic field, behaving niz ^ \h — hc\^ with 
he = 0.604, l3 = 0.98(1) for a = 0.3 and D = 5. Our re- 
sult is consistent with the theoretical prediction /3 = 1.0 in 2D 





FIG. 4: (Color online) The ground state phase diagram of the spin- 
1/2 HAF system on the DHC lattice in a-h plane, where three phases 
(dimer, semiclassical Neel and polarized) are identified, and the 
collinear Neel phase marked by the star line inhabits exactly on the 
a axis (h = 0). The fitting curves reveal the critical behavior in the 
region a e [0.45,0.51]. 
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FIG. 5: (Color online) The log-log plot of spin gap A versus a in the 
critical region a £ [0.45, 0.51], where the linear fit gives the critical 
exponent u ~ 0.69(2). 



and higher quantum spin systems*'. Nonetheless, it is in sharp 
contrast to the cases of gapped one-dimensional Heisenberg 
spin systems where it is observed a square root dependence 
of niz ^ \h — hc\^^'^ that characterizes the commensurate- 
incommensurate phase transition^*. For other a < ac, such 
an almost linear dependence was also noted. The suscepti- 
bility, X — dniz/dh, as a function of the magnetic field h 
is shown in the upper inset of Fig. [3] (a). Two discontinu- 
ous points in x versus h are seen, namely one at the point he 
where the spin gap closes, and the other at the saturation field. 
We have also explored the transverse component of sublattice 
magnetization per site, ttIx^a — j^J^ieA^^f)' against the 
magnetic field h, as given in the lower inset of Fig. |3ja). For 
a < ac, rnx,A vanishes for h < he, while it increases sharply 
when h exceeds he, and after reaching a round peak it de- 
clines steeply and vanishes at the saturation field. This obser- 
vation shows that, with increasing the field, there exists a tran- 



FIG. 6: (Color online) The staggered magnetization per site ml as 
a function of the staggered magnetic field hs for h = 0. The inset 
is 1/Xs against a, where the fitting curve is obtained at a small field 
hs = O.OOIJ. 



sition from the disordered dimer phase to a canted Neel phase 
(spin flop phase)'*'^'^^ that is characterized by nonzero values 
of both rriz and rrix^A, where the spins align antiferromagneti- 
cally within the xy plane and develop a uniform z component 
along the field, thus canting out of the plane. For a > ac, 
the magnetization curves behave differently from those with 
a < ac, as presented in Fig. [3](b) for a = 0.8 as an example. 
With increasing the magnetic field, the magnetization m^ in- 
creases monotonously till the saturation, while the transverse 
component of magnetization TOj.^ first increases slowly, then 
drops sharply and vanishes eventually at the saturation field. 
The susceptibility x increases slowly with increasing the field, 
and then decreases steeply to zero at the saturation field. 

By summarizing the above observations, a ground-state 
phase diagram of the system in the a — h plane can be drawn, 
as presented in Fig. HI where the phase boundaries are deter- 
mined by the transition points in Figs. |2]and[3] One may see 
that there are three phases, namely, the dimer phase, the semi- 
classical Neel phase (including the canted and collinear Neel 
states) and the polarized phase. At a = ac, there is a QPT 
from the disordered dimer phase to the collinear Neel phase. 
Note that the lower phase boundary between the dimer and 
canted Neel phases is determined by observing the spin gap A 
that is obtained by calculating the width of zero magnetization 
plateau for various a as presented in Fig. [3] (a). The critical 
behavior of the spin gap in the dimer phase in the vicinity of 
ac is fitted by the least squares method with A ^ (a — acY, 
where ac — 0.54, as shown in Fig. |4] The critical exponent is 
found to be 1/ ~ 0.69(2) by the linear fit, as shown in Fig. |5] 
It is close to the standard 0(3) value of 0.7112(5) (e.g. Refs. 
Ial9ll27l) . In comparison to the result ac = 0.27 of the nonlin- 
ear a model method^^ and the variational result a^var = 0.5 
obtained through Eq. (O, the present TRG result is closer to 
ac = 0.576 and v = 0.707 of QMC calculations' I As ex- 
pected, owing to its lower coordinates, the disordered region 
on the DHC lattice is wider than that on a square lattice where 
ac = 0.397i2. 

In the presence of a staggered magnetic field hg, the stag- 



gered magnetization per site, ml — {^/NjJ^i^iiSi), as a 
function of hg for different a is presented in Fig. |6] where 
h — 0. With increasing hg, ml increases monotonously 
from zero for a < a^, while for a > ac, ml starts from 
a nonzero value, implying again that the system has a spon- 
taneous AFLRO for larger a. The inverse of staggered sus- 
ceptibility Xs — dml/dhs as a function of a is presented 
in the inset of Fig. |6]for h — and a very small hg. The 
critical behavior of the zero-field staggered susceptibility Xs 
is expected to diverge as \a — ad^^ near the critical point 
ac, where the nonlinear curve fitting gives Uc — 0.542(1) 
and the log-log plot in the inset shows the critical exponent 
7 ~ 1.363(8), which agree with Uc obtained through the spin 
gap and the 0(3) value of 7 = 1.373(3)^^, respectively. Ac- 
cording to the obtained critical exponents (3 and 7, it is seen 
that the quantum criticality of the present system falls into the 
0(3) universality class. Similar calculations indicate that the 
transition between the dimer and canted Neel phases have the 
same critical exponents j3 and 7 as that of the dimer-collinear 
Neel phase transition at h = 0. Hence, they also belong to the 
same 0(3) universality class. 



V. RELATION TO THE 2D BOSON HUBBARD MODEL 

Finally, we would like to mention briefly that the present 
spin- 1/2 system has a link with the 2D boson Hubbard model. 
By performing the hard-core boson mapping^'^ with aj -^ S^, 
Ui — >■ S^, rii — ^ Sf + 1/2, and making a rotation of spins by 
■K along the z axis on one sublattice, one can obtain, from Eq. 
(1) with hs = 0, the Hamiltonian of hard-core bosons 

Hb = Y^ [-t{a\aj + a]a{) + U{n, - l/2){nj - 1/2)] 

+ Yl [-^'(«l«i + «!«») + ^'("« ~ l/2)("j - 1/2)] 
- Ai^(n,-l/2), (6) 

i 

where t = U = J, t' = U' = J' and /i = /i. It is the boson 
Hubbard model on the DHC lattice, whose behavior can be 
well understood in accordance with the aforementioned cor- 



responding Heisenberg spin system. With the above mapping, 
owing to the nonvanishing sublattice magnetization |toj. ^(s) | 
in x-y plane (see the inset of Fig. |3]l, the canted Neel phase 
in the spin system corresponds to the boson superfluid phase 
with off diagonal long-range order (ODLRO)^''^" in the boson 
Hubbard model for i'/t = U'/U > 0.54 on the DHC lattice. 
For t' /t = U'/U < 0.54, the disordered dimer phase in the 
spin system is mapped onto a liquid state with neither ODLRO 
nor DLRO. Therefore, the boson Hubbard model on the DHC 
lattice has a QPT from the boson liquid to a superfluid. The 
spin polarized phase becomes a Mott insulator phase with one 
boson occupying each site in the mapped boson system, and 
hence, there also exists a superfluid-Mott insulator transition. 
VI. CONCLUSION 

In conclusion, the spin- 1/2 HAF model on the DHC lattice 
is studied by means of the combined iTEBD and TRG algo- 
rithm, where the ground-state phase diagram is obtained. It 
is uncovered that there is a second-order QPT from the dis- 
ordered dimer phase to the ordered collinear Neel phase at 
ac — 0.54, where the critical exponents v ~ 0.69(2) and 
7 ~ 1.363(8) are determined. This QPT belongs to the stan- 
dard 0(3) universality class. In addition, through a boson 
mapping, the present HAF system has a link with the boson 
Hubbard model on the DHC lattice. The properties of the lat- 
ter boson Hubbard model can thus be understood in terms of 
the present study. We expect that our findings is not only use- 
ful for understanding experimental observations of the anti- 
ferromagnets with DHC lattices, but also is helpful for the 
corresponding 2D boson Hubbard model. 
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